www.gusucode.com > C-C方法计算时间延迟和嵌入维数计算Lyapunov指数计算关联维数混沌时间序列预测 > code/Chaos Toolbox Ver.2.0/Main/MainPre_by_Lya_n.m

    %skyhawk
clear all;
A=load('bk1.txt');
whl=A(:,4);
[whsl,lll]=size(whl);  

m=6

P=26;  %平均循环周期

N=80;  %最多校验步数

for i=1:whsl
    whlsj(i)=whl(i);
end

[lmd_m,idx,min_d,idx1,min_d1]=lyapunov(m,whlsj,whsl,P);%求lyapunov指数
lmd_m
% t_m=fix(1/lmd_m)       %最大预测步数
t_m=5

for i=(whsl-N+2-t_m):(whsl-N+1)         %预测后t_m步
    [lmd_m,idx,min_d,idx1,min_d1]=lyapunov(m,whlsj,whsl-N+2-t_m-1,P);  
    [yc,z(i)]=pre_by_lya(m,lmd_m,whlsj,i-1,idx,min_d);%预测第i+1个点  
    whlsj(i)=yc;  %将第i步预测值作为完好率数据的第i个值进行下一步预测
end

y(whsl-N+1)=yc;
fch=(y(whsl-N+1)-whl(whsl-N+1))*(y(whsl-N+1)-whl(whsl-N+1));
shuliang=1;

for i=(whsl-N+2):(whsl)
    whlsj(i-t_m)=whl(i-t_m);  %换为实际值
    [lmd_m,idx,min_d,idx1,min_d1]=lyapunov(m,whlsj,i-t_m,P);  
    [y(i),z(i)]=pre_by_lya(m,lmd_m,whlsj,i-1,idx,min_d);%预测第i+1个点  
    whlsj(i)=y(i);           %换为预测值
    
    fch=fch+(y(i)-whl(i))*(y(i)-whl(i));
    shuliang=shuliang+1
end

fch=sqrt(fch)/shuliang

yyy=[whl,y'];
save('kjycqush.txt','yyy','-ASCII');

kk=1:whsl;
plot(kk,whl,'b',kk,y,'r')